import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

def set_chinese_font():
    """设置中文字体"""
    plt.rcParams['font.sans-serif'] = ['SimHei', 'Microsoft YaHei', 'DejaVu Sans']
    plt.rcParams['axes.unicode_minus'] = False  # 解决负号显示问题

set_chinese_font()

def mathematical_modeling():
    # 示例1: 人口增长模型
    def population_growth_model():
        # 生成模拟数据
        years = np.arange(2000, 2021)
        # 使用指数增长模型: P(t) = P0 * e^(rt)
        base_population = 1000
        growth_rate = 0.02
        population = base_population * np.exp(growth_rate * (years - 2000))
        
        # 添加一些随机噪声
        np.random.seed(42)
        population_noisy = population + np.random.normal(0, 50, len(years))
        
        # 定义指数模型函数
        def exp_model(t, P0, r):
            return P0 * np.exp(r * (t - 2000))
        
        # 拟合模型
        popt, pcov = curve_fit(exp_model, years, population_noisy, p0=[1000, 0.02])
        
        # 预测未来人口
        future_years = np.arange(2000, 2031)
        predicted_population = exp_model(future_years, *popt)
        
        # 可视化
        plt.figure(figsize=(10, 6))
        plt.scatter(years, population_noisy, color='blue', label='实际数据')
        plt.plot(future_years, predicted_population, 'r-', label=f'拟合模型: P(t) = {popt[0]:.1f} * e^({popt[1]:.3f}t)')
        plt.axvline(x=2020, color='gray', linestyle='--', alpha=0.7, label='当前年份')
        plt.xlabel('年份')
        plt.ylabel('人口')
        plt.title('人口增长模型')
        plt.legend()
        plt.grid(True)
        plt.show()
        
        print(f"人口增长模型参数:")
        print(f"  初始人口 P0 = {popt[0]:.2f}")
        print(f"  增长率 r = {popt[1]:.4f}")
        print(f"  2030年预测人口: {exp_model(2030, *popt):.0f}")
    
    # 示例2: 抛物线运动模型
    def projectile_motion():
        # 初始条件
        v0 = 50  # 初始速度 m/s
        angle = 45  # 发射角度 度
        g = 9.8  # 重力加速度 m/s²
        
        # 转换为弧度
        theta = np.radians(angle)
        
        # 计算飞行时间
        t_total = 2 * v0 * np.sin(theta) / g
        
        # 生成时间点
        t = np.linspace(0, t_total, 100)
        
        # 计算位置
        x = v0 * np.cos(theta) * t
        y = v0 * np.sin(theta) * t - 0.5 * g * t**2
        
        # 可视化
        plt.figure(figsize=(10, 6))
        plt.plot(x, y, 'b-', linewidth=2)
        plt.xlabel('水平距离 (m)')
        plt.ylabel('高度 (m)')
        plt.title('抛物线运动模型')
        plt.grid(True)
        
        # 标记关键点
        max_height = (v0 * np.sin(theta))**2 / (2 * g)
        max_range = v0**2 * np.sin(2*theta) / g
        
        plt.axhline(y=max_height, color='r', linestyle='--', alpha=0.7, label=f'最大高度: {max_height:.1f}m')
        plt.axvline(x=max_range, color='g', linestyle='--', alpha=0.7, label=f'最大射程: {max_range:.1f}m')
        plt.legend()
        plt.show()
        
        print(f"\n抛物线运动模型:")
        print(f"  最大高度: {max_height:.2f} m")
        print(f"  最大射程: {max_range:.2f} m")
        print(f"  飞行时间: {t_total:.2f} s")
    
    population_growth_model()
    projectile_motion()

mathematical_modeling()